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OPTIMIZATION OF SYSTEMS WITH UNCERTAINTY: INITIAL DEVELOPMENTS 
FOR PERFORMANCE, ROBUSTNESS AND RELIABILITY BASED DESIGNS* 

LUIS G. CRESPO* 


Abstract. This paper presents a study on the optimization of systems with structured uncertainties, 
whose inputs and outputs can be exhaustively described in the probabilistic sense. By propagating the 
uncertainty from the input to the output in the space of the probability density functions and the moments, 
optimization problems that pursue performance, robustness and reliability based designs are studied. By 
specifying the desired outputs in terms of desired probability density functions and then in terms of mean- 
ingful probabilistic indices, we settle a computationally viable framework for solving practical optimization 
problems. Applications to static optimization and stability control are used to illustrate the relevance of 
incorporating uncertainty in the early stages of the design. Several examples that admit a full probabilistic 
description of the output in terms of the design variables and the uncertain inputs are used to elucidate 
the main features of the generic problem and its solution. Extensions to problems that do not admit closed 
form solutions are also evaluated. Concrete evidence of the importance of using a consistent probabilistic 
formulation of the optimization problem and a meaningful probabilistic description of its solution is provided 
in the examples. In the stability control problem the analysis shows that standard deterministic approaches 
lead to designs with high probability of running into instability. The implementation of such designs can 
indeed have catastrophic consequences. 

Key words, structured uncertainty, optimization, robustness, reliability, performance, control 
Subject classification. Applied and Numerical Mathematics 

1. Introduction. Most of the engineering analysis tools rely on a precise mathematical description 
of the physical problem. Once a model is built, deterministic procedures are usually applied to study the 
behavior of the physical system. Conjectures are made and decisions are taken based on the results of this 
process. Nevertheless, in most of the practical cases we can not specify precisely the model (unstructured 
uncertainty) or the value of its parameters (structured uncertainty). The assumption of deterministic values 
for the uncertainty might lead to faulty models whose analysis might not be in conformity with reality. 

This paper focuses on problems with structured uncertainties. Parameter uncertainties are typically 
specified in terms of interval analysis, memebership functions and probability density functions (PDFs) [10]. 
The main characteristic of interval analysis is that variables are represented by lower and upper bounds. 
This approach has been applied in a variety of fields, including robust controls [1, 2, 13], structural analysis 
[16] and trajectory planning [11]- This is the least accurate method for uncertainty modeling. When the 
parameter uncertainties are characterized by membership functions, fuzzy logic is the basis for assessing the 
uncertainties in the system’s output [4]. This description provides an intermediate level of detail. 

Uncertainty-based design methods using PDFs are referred to as probabilistic methods. Such methods 
provide the best description of the uncertain parameters by treating them as random variables [8]. Monte 
Carlo Simulation, Importance Sampling, Latin Hypercube Sampling and Generalized Cell Mapping are 
numerical methods commonly used to estimate PDFs [7]. Furthermore, stochastic differential equations are 
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used to propagate uncertainty in initial conditions, boundary conditions, and excitations. Stochastic controls 
[3] and Polynomial Chaos [18] typify this practice. 

Optimization of systems with uncertainty presents further challenges. In addition to the description 
and propagation of the uncertainty, the output is now parametrized by the design variables (deterministic 
quantities). In other words, for every possible design the uncertainty must now be propagated from the 
input to the output. This obviously makes the problem much more resource-intensive. In addition, opti- 
mization problems with different objectives such as performance, robustness and reliability must be properly 
formulated and solved in order to provide meaningful solutions. 

The high degree of difficulty in obtaining refined probabilistic descriptions of inputs and outputs has 
obstructed the application of existing mathematical tools. Considerable effort has been devoted to the 
development and implementation of viable numerical strategies. Sampling-based methods constitute the 
core of such approaches. By fitting all available data with some smooth surface, i.e. response surface, a 
model of the system is built. Then, standard optimization schemes are applied. Response surface methods 
[17] and Taguchi parameter design methods [9] are commonly used in practice. 

This paper studies the optimization of systems with uncertain parameters by using exact expressions 
for the PDFs and the moments of the output. This treatment allows to attain an exhaustive probabilistic 
description of the output/response in terms of the design variable and the statistical properties of the input. 
In this framework, optimal performance, robustness and reliability based designs can be found without any 
tergiversation. 

The content is organized as follows. Section 2 presents some basic mathematical tools required for further 
developments. Definitions of engineering design criteria and their corresponding metrics are introduced in 
Section 3. Section 4 presents several applications to static optimization and stability control. The solution to 
selected examples can be used not only to evaluate the accuracy of numerical methods but more importantly 
to understand the very primary nature of the problem. By stating some remarks at the end of each section 
we intend to build up some insight and understanding of the qualitative features of the problem at hand, its 
proper formulation, and its solution. Finally, some conclusions are stated in Section 5. 

2. Background. Let a be a random variable with cumulative probability distribution F A (a) and prob- 
ability density function f A (a). If the functional relation between a and y is given by y = Cl{a), the PDF of 
y is given by: 

fy(y) = ^P{ft(A) <y} = i p ( A ® (2.1) 

where r(y) represents the values of a such that ft (.4) = y and refers to the inequality symbol that 
corresponds to ft(-4) < y. This event might be formed by the combination of several mutually exclusive 
events. Before differentiating, the event A ® r(y) must be fully written in terms of F A (a). 

When y is a function of n random variables, an n-dimensional integral must be calculated to find fy(y ). 
In particular, if y — ft(ai,... ,a„) and f Au ... ,A n (ai , . . . ,a n ) is the joint probability density function of the 
random variables A \, . . . , .4 n , the PDF of y is given by: 

f y { y ) = ~r f ... f f A\ 4„(ai,... ,a n )dai ...da n , (2.2) 

The same results can be obtained via the characteristic function [5]. 

3. Definitions and Metrics. Throughout this manuscript, we denote with x the design variable, a the 
uncertain input parameter and y = ft(x, a) the system output. For our purposes the uncertain parameters are 
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Fig. 3.1. Propagation of probability density functions through the operator Q(x, a) 


random variables distributed according to prescribed density functions. See Figure 3.1. By propagating the 
PDF of the uncertain parameters through the system we can attain an exhaustive probabilistic description 
of the output via its PDF. Once the effects of the design variable on the output are fully understood in the 
probabilistic sense, the optimization can be carried out. In this context, the optimization process intends to 
shape the PDF of the output according to prescribed objectives. 

Designs with performance, robustness and reliability specifications are needed in practice. Let’s first 
define conceptually what such designs intend to achieve [19]. A Performance-based design intends to minimize 
the expected value of the cost while taking into consideration the entire ensemble of values that the uncertain 
parameters might take. A Robustness-based design intends to be insensitive to changes in the uncertain input 
variables. In other words, a robust design intends to respond in the same fashion for the whole set of values the 
uncertain parameter might take. A Reliability-based design seeks to minimize the probability of occurrence 
of a catastrophic event. Clearly, these criteria define different optimization problems. Let’s now state the 
following basic questions: 

• How should the optimal PDF of the output be for each of these designs? 

• How do we properly formulate the corresponding optimization problems? 

• Do these output specifications lead to conflicting designs? 

This paper intends to elucidate these matters by solving optimization problems that have explicit math- 
ematical solutions. Notice however, that finding closed form solutions for the output is not possible in most 
of the practical problems. The high degree of difficulty of obtaining PDFs for the output forces us to down- 
grade the quality of its probabilistic description. By taking into account the true features of the problem at 
hand along with its practical impediments, consistent simplifications must be made. 

With this spirit we define a set of indices whose minimization leads to the desired optimal designs. The 
indices for performance, robustness and reliability based designs to be used are Var [j/] and P{y > yum } 
respectively. The reader must notice that these metrics capture the essence of the desired output in viable 

probabilistic terms. 

We will focus on problems whose desired output does not combine performance, robustness and reliability 
requirements. The reasons for this will be apparent later on. Nevertheless, it is important to highlight that 
a practical problem usually has a combination of performance, robustness and reliability requirements. 
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4. Applications. In this section we study the propagation of the PDF of the uncertain input(s) to the 
output. The optimal performance, robustness, and reliability based design problems are studied by solving 
selected examples in static systems and stability control. Remarks of generic scope give closure to each 
application. 

4.1. Static Optimization. The problem formulation to all the examples in this section is as follows. 
Assume a and b are input random variables distributed according to f A (a) and f 11 (b), x is a design variable 
and y is the output of the system. If the input(s)-output relation is given by y = fi (x,a,b), how does fy(y) 
vary with respect to x and what are the optimal designs? 

The optimal solution to the problem in which the uncertain parameters take the value of their mean with 
probability one is called the deterministic optimal solution. The following notation is used in the examples. 
For performance-based design x* e = argmin{ £[,</]} and y* e = inf, {£%]}, for robustness-based design x’ T0 = 
argmin{Var[t/]} and y’ 0 = inf x {Var[j/]} and for reliability based design x* e = argmin{P{.y > y lim }} and 
y re = inf x {P {y > yiim }}• 

Examples that admit closed form expressions for f Y (y) are presented first. Notice that the resulting 
expressions apply to any distribution of the input. Exact expressions for the first two moments of the output 
are also derived. In addition, numerical schemes that allow us to extend the existing mathematical tools to 
systems that do not admit exact solutions are explored. 

4.1.1. Ideal Robust Design. Consider the system y = fi(a, x) = a 2 x 2 4- ax. 

The deterministic solution of this problem leads to argminjy} = -l/( 2a) and inf x {?/} = -1/4. Notice 
that when a moves through zero, argminjj/} moves towards infinity, where it changes its sign. If the PDF 
of the input is given by f A (a ), the system output is distributed according to: 

fr(y) = Csign(x){l/{kx))[f A ((-l 4 k)/2x) 4 f A {{- 1 - k)/2x)\, (4.1) 

where C is a normalization constant and k = y/1 4 4x . Notice that for real values of y , a bound on the 
range of fy{y) is imposed. This bound is in agreement with the inf x {y} in the deterministic problem. The 
spectrum of PDFs for all x € [—1.8, 1.8] when a is a Gaussian random variable with mean 0.5 and standard 
deviation 0.6, i.e. a iV(0.5,0.6), is shown in Figure 4.1. Moments of arbitrary order of y and probabilities 
of failure can be readily found by integrations that involve Equation (4.1). 

This approach by itself is unable to provide the explicit dependence of the moments of y on the parameters 
of /a(g) and x. On the other hand, the moment generating function for the PDF of a can be systematically 
used to find this information. In this work we use both approaches to extract optimal solutions for the three 
classes of problems we are interested in. 

If a *4 A r (m, 5), the mean and variance of y are given by: 

E[y] — (m 2 4 s 2 )x 2 4 mx (4.2) 

Var[y] = (4m 2 s 2 4 2s 4 )x 4 4 4ms 2 :r 3 4 s 2 x 2 . (4.3) 

The optimal solutions of the three designs are as follows: 

> Performance : x* e = -m/[2(m 2 4 s 2 )] and y* e = -m 2 /[4(m 2 4 s 2 )] 

> Robustness : x* Q = 0 and y* 0 = 0 

> Reliability: x* e = -0.2824 and y* e - 0.22. 

While the first two solutions are found using explicit expressions for the moments of y , the last one is found 
via fy(y) using yu m = —0.01. As expected, results based on Equations (4.2) can also be obtained using 
Equation (4.1). 


4 




Fig. 4.1. Spectrum of PDFs of the output for a Gaussian input random variable 



Fig. 4.2. Optimal deterministic (dashed line), performance (solid line), robustness (dotted line) and reliability (dash- dotted 
line) PDFs of the output for a Gaussian input. A vertical dotted line at yn m is also shown . 


The corresponding optimal PDFs are shown in Figure 4.2. Notice that the total elimination of the 
variance is achieved by multiplying the uncertain parameter by zero. This is possible due to the particular 
structure of D(x,a). 

The optimal robust design leads to a delta function at zero for the PDF of y. This is the ideal case 
because the system output is completely insensitive to variations in the uncertain parameter. Notice also, 
that the solution implies that f Y (y) = 0 for y < -1.4. This fact is in agreement with the corresponding 
deterministic problem. 

4.1.2. Bounds in the PDFs of the Output. Consider the system y = Q{x,a) = x 2 + ax + a 2 . 

The deterministic solution to this problem leads to argmin{i/} = — a/2 and inf x {y} = 3a 2 /4. In contrast 
to the previous example, the design variable is unable to cancel the effect of uncertainties in a. Hence, an 
optimal robust design would reduce the spread of fy(y) about its mean value without fully eliminating it. 

If a is distributed according to /^(a), the PDF of the output is given by: 

f Y [y) = C(l/k)[f A ((-x + k)/ 2) + f A ((-x - k)/2)\, (4.4) 

where C is a normalization constant and h = >/4y — 3x 2 . The optimal PDFs for a — > 7V(0.6,0.2) and 
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FlG. 4.3. Optimal deterministic (dashed line), performance (solid line), robustness (dotted line) and reliability (dash-dotted 
line) PDFs of the output for a Gaussian input. A vertical dotted line at yu m i s also shown. 



Fig. 4.4. Optimal deterministic (dashed line), performance (solid line), robustness (dotted line) and reliability (dash-dotted 
line) PDFs of the output for a Beta input. A vertical dotted line at yum is also shown 

a -4 £(1.2, 1.2) are shown in Figures 4.3 and 4.4. Here £(1.2, 1.2) denotes a beta distribution that is 
symmetric about a = 0.5 and has a bounded support given by a E [0,1]. The optimal reliability-based 
designs were searched within the interval x € [-2.5, 1.5]. 

If a — ► N(m, s), the first two moments of y are given by: 

E[y\ = x 2 4- mx 4 m 2 4 s 2 (4.5) 

Var[y] = s 2 x 2 4 4ms 2 x 4 4 m 2 s 2 4 2s 4 . (4.6) 

Optimization leads to the following results: 

> Performance : x* e = — m/2 and y* e = s 2 4 3m 2 /4 

> Robustness: x* 0 = —2m and y* 0 = 2s 4 

o Reliability: x* Te — —0.409 and y* e = 0.141 

where we assumed yu m — 0.5. Because y is real, a bound on the domain of fy(y) is imposed regardless of 
the shape of /^(a). In contrast to the previous example, this bound, i.e. y > 3x 2 /4, depends on the design 
variable. It is interesting to notice that the shapes of fy(y) for x* c , x* ro and x* re are similar to the ones found 
in the previous example. This observation also applies to the examples to come. Explicit expressions for the 
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2 . 5 , 



Fig. 4.5. Spectrum of PDFs of the output for a Gaussian input random variable 


moments of y can also be derived for other forms of f A (a)- For instance, if a — > B(p,q) the optimals are 
given by x* e = -p/2(p + q) and x* ro = -2(1 + p)/(2 + p + q)- 

4.1.3. Multiple Extrema. Consider the system y = fl(x,a) = x 4 /4 + (1 - a)x 3 /3 - ax 2 / 2 4- a 2 . 

The deterministic solution leads to three extrema for y. They occur atx = -l,x = 0 and x = a. We will 
refer to these points as Xj, X2 and X3, respectively. After some manipulations we find. 

y(xi) = a 2 -a/6 -1/12 
y(x 2 ) = a 2 

y(x 3 ) = a 2 (l — a/6 — a 2 /12) 


argmin{y} = < 


X\ if a € [-0.5, 1] 
X 2 if a E [—2, —0.5] 


^ £3 otherwise. 


(4.8) 


If a is distributed according to /.4(a), the PDF of the output is given by: 

Ml/) = C(l/k)[f A ((g + k)/ 2) + f A ((g - *0/2)], ( 4 - 9 ) 

where C is a normalization constant, y = x 2 / 2 + x 3 /3 and k = y/g 2 + 4(y - x 4 /4 - x 3 /3). The spectrum 
of PDFs for all x G [-2,2] with a -y 7V(1,0.4) is shown in figure 4.5. The corresponding optimal PDFs are 

shown in Figure 4.6. 

In this case, the bound on f Y (y ) is given by the positiveness of k. The high-order polynomial dependence 
of the bound on the design parameter is a consequence of the non-quadratic structure of fl(x, a). Notice that 
the non-trivial expression for this bound appears naturally in the process. The existence of several extrema 
in the deterministic problem might lead to non-unique optimal solutions in the probabilistic problems. 

If a -4 N(m, s), the first two moments of y are given by: 

E[y] = x 4 /4 + x 3 /3 — mg + m 2 + s 2 ( 4 - 10 ) 

Var[y] = s 2 (x 6 / 9 + x B /3 + x 4 /4 - 4mx 3 /3 - 2mx 2 + 2s 2 + 4m 2 ). (4.11) 


Optimization leads to the following results: 
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Fig. 4.6. Optimal deterministic (dashed line), performance (solid line), robustness (dotted line) and reliability ( dash-dotted 
line) PDFs of the output for a Gaussian input. A vertical dotted line at y lim is also shown. 



Fig. 4.7. Parametrical curves for the first two moments of y for Var[a] - 0.01 and a varying E[a] 


> Performance: x* pe = argmi n{y} and y* e = y{x* pe ) + s 2 taking a = m 

> Robustness: x* 0 = xj and y* 0 = Varfzj] 

> Reliability: x* re = 1.248 and y* e = 0.2265, 

where xj must be selected out of x u x 2 and the real roots of 2x 3 + 3x 2 - 12 m = 0. The explicit expression 
for Xj is omitted due to space limitations. For the reliability problem yu m = 1.2 was used. 

From these results, several observations can be made. The offset between the deterministic solution and 
the optimal performance based design is just a constant. In addition, while the variance depends upon both 
m and s, the location of its extrema depends on m exclusively. 

By eliminating the design variable from Equations (4.10) and (4.11), we generate parametric curves that 
relate the first two moments of the output. Once the mean and the variance of a are set, a single parametric 
curve applies. For a fixed PDF of the input parameter, there is a one-to-one correspondence between one 
design and one point on this curve, i.e. changing the value of x moves the operating point on the curve. 
In this fashion, we can study exactly and simultaneously the sensitivity of both the mean and the variance 
of the output to changes in the design variable. Furthermore, we can evaluate the effects of formulating 
optimization problems with just one of the indices. 
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Fig. 4.8. Zoom of the previous figure 



Fig. 4.9. Parametrical curve for the first two moments of y for E[a] — 0.5 and a varying Var[a] 


Figure 4.7 shows the parametric curves for s — 0.1 and m varying from —0.25 to 0.35 by increments of 
0.1 units. The dotted line corresponds to m = 0.35. We can see that for m < -0.2, the minimization of 
incidentally increases the variance. This behavior is clearly undesirable in practical applications. This 
simple case provides evidence of the need for formulations that intend to reduce the variance of the output 
regardless of the type of design we are interested in. 

A zoom of the region of minimum variance is shown in Figure 4.8. From this plot we can see that the 

solution to the minimum variance problem might not be unique. 

In Figure 4.9, a set of parametric curves in which the mean of a is kept constant while its variance is 
increased are shown. In this case, m = 0.5 and s is increased uniformly from 0.1 with increments of 0.1 
units. The dotted line is the curve with the maximum value of s. The reader must notice that both, E[y] 
and Var[y] have more than one extremum. From this figure we can see that an optimal robust design and 
an optimal performance design might differ considerably. 

Figure 4.10 shows the parametric curves for a varying input variance s, when the mean is kept constant 
at m = 1. From this figure, we can see that two different designs lead to the same optimal performance 
while the corresponding variances are quite different, specifically for large values of s. In this case the higher 
order moment can be used to discriminate between these two designs. 
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Fig. 4.10. Parametrical curve for the first two moments of y for E[a] ~ 1 and a varying Vai\a] 



Fig. 4.11. Spectrum of PDFs of the output for a Gaussian input random variable 

4. 1.4. Non-polynomial Uncertainties. Consider the system y = ft(z,a) = y/x 2 + ax + a 2 . 

The argument of the square root is the function used in Example (4.1.2). The deterministic solution is 
given by argminji/} = —a/2 and inf z {2/} = |a|-y/3/4. For a given /4(a), the system output is distributed 
according to: 


fy(y) = C(y/k)[f A ((-x + k)/ 2) + f A ((-x - k)/ 2)], (4.12) 

where C is a normalization constant and k = ^jAy 2 -3x 2 . In this case the bound on E[y] is given by 
V > \ x \y/3p- The spectrum of PDFs for all a; € [-2.5, 2.5] is shown in Figure 4.11. 

Figure 4.12 shows the PDFs that correspond to the optimal designs. In the calculation we assume 
Viim = 1- Recall that this approach provides exact results. 

The fact that fi(x,a) is not polynomial on a prevents us from obtaining exact and explicit expressions 
for the moments of y in terms of the parameters of f A (a) and x. Given the relevance of this information 
to the optimization problem and the non-polynomial structure of most of systems, two schemes to calculate 
the moments of y as a function of x axe used next. The first one is via polynomial approximations of n(ar, a) 
and the second one is by approximating the expressions for the moments of y. 

Approximations of fi(x,a): Assume that the function f l(x,a) approximates well the function fl(x,a) in 
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Fig. 4.12. Optimal deterministic (dashed line), performance (solid line), robustness (dotted line) and reliability (dash- 
dotted line) PDFs of the output for a Gaussian input A vertical dotted line at y lim is also shown. 

the interval of interest. More precisely, 

n 

H (x,a) ~ Q(z,a) = ^ C{(x)a l for a € [ai,a 2 ], (4.13) 

t=0 

where n can be arbitrary large. This approximation can be done using a Taylor expansion about E[a\ or 
using an interpolating polynomial that satisfies Q(x, Oj) = Q(x, a*) for i = 1 . . . n. Once the approximation is 
built, the expected value operator can be applied directly. Notice that nothing prevents us from using high 
order polynomials. This approach provides an approximation of the explicit dependence of the moments of 
the output on the parameters of /^(a) and x. 

The quality of the Taylor expansion might be severely diminished with large deviations from the point 
of expansion. This makes it suitable for cases where /^(a) is highly concentrated about its mean, i.e. 
Var[a] <£ 1. On the other hand, the range in which an interpolating function approximates well fi(x,a) is 
given by the values a* used to build it. This feature suites well problems in which the PDFs of the input are 
dispersed. 

Notice that for a Gaussian /^(o), the approximation of fi(x,a) must be good in the range a e [m - 
3s, m + 3s] for all x. Due to the exponential decay of the PDF, large offsets outside this range will not affect 
the resulting expressions for the moments. Considerations of this type can be used to restrict the sampling 
space according to the relevant intervals of the support of the PDF of the inputs. 

Approximation of the moments: This method is mounted on the following approximation: 

/ oo roo 

y m f Y {y)dy= / n(x,a) m f A (a)da ~ ^ D(x, (4.14) 

-oo J—oo i— o 

where the summation is a Riemann sum finely partitioned on the interval in which f A (A) is non-zero. In 
practice, this summation can be evaluated with any numerical scheme. While derivatives of Q(x,a) are not 
required, a fine partition of it is needed for accurate results. If non-uniform partitions are used, they must 
be particularly fine where a has a higher chance to occur and D(x,o) reaches its extrema. Notice that higher 
chance of occurrence does not necessarily happen in the vicinity of £[o]. As in the interpolating polynomial 
scheme, this method has the advantage of not concentrating the entire approximation about a single point. 
Although this method does not provide the explicit dependence of the moments on x, its computational 
simplicity makes it practical. 
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Fig. 4.13. Exact and approximate curves for the first two order moments of y as a function of x: Exact (solid line), 
Taylorfdashed line) and moments approximations for uniform (dash-dotted line ) and non-uniform partitions (dotted line) 


These approaches are applied to the particular problem of interest. Figure 4.13 shows the exact and 
approximate results for E[y] and Var[y] as a function of x when a -4 7V(0.7,0.3). Optimal designs can be 
extracted from this information. While in this case, all numerical approaches lead to satisfactory results, 
the method based on the approximation of the moments performs the best. 

For all these approximations, notice that the offset between the true value of the moments and its 
approximation is propagated and accumulated through the calculation of higher order moments. 

4.1.5. Multiple Uncertain Parameters. Consider the system y = Q(x, a, 6) = x 2 4 ax 4- 6 2 , where 
a -4 N(m,s) and b -4 jV(/x, a) are independent random variables. 

Due to the high degree of difficulty in deriving exact and approximate PDFs for the output for multiple 
uncertain inputs, we rather focus on finding exact expressions for its moments. This practice not only allows 
us to solve for the performance and robust optimal designs but also is easily extendable to problems with 
several independent inputs. By applying the expected value operator to ft(x , a,5) we find: 

E[y] = x 1 4 mx + fj, 2 + a 2 
Var[j/] = s 2 x 2 4 Aft 2 a 2 4 2a 4 

The optimal designs are given by: 

> Performance : x* e = -m/2 and y* e = fi 2 4 o 2 — m 2 /4 

> Robustness : x* Q = 0 and y* Q = Var [b 2 ] 

As in Example 4.1.1, the optimal robust design is independent of a. 

Let's now consider y = Q(x, a, b) = x 2 4 xa 2 4(1- x)b 2 . In this case the design variable faces a trade-off 
between both input variables. The corresponding expressions for the moments are: 

E[y\ - x 2 4 x (m 2 4 s 2 - /i 2 - a 2 ) 4 fi 2 4 a 2 (4.17) 

Var [y] = 2 rjx 2 - 4<r 2 £x 4 2cr 2 f , (4.18) 

where i)e(s 4 4<j 4 4 2m 2 s 2 4 2/i 2 a 2 ) and £ = o 2 4 2 ft 2 . Optimization leads to: 

t> Performance : x* pe = (ft 2 4 o 2 - m 2 - s 2 )/ 2 and y* € =/i 2 4 ( 7 2 - (m 2 4 s 2 - fi 2 - a 2 ) 2 / 4 

> Robustness : x* 0 = a 2 ^/rj and y* D = 2s 2 a 2 (2m 2 4 s 2 )£/rj. 


(4.15) 

(4.16) 
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Fig. 4.14. Parametric curves for the first two moments of y as the variance of one of the parameters vanes 


In contrast to the previous problem, the robust solution is now not trivial. The corresponding parametric 
curves for m = 0.5, n = -0.6, a = 0.5 and a varying s are presented in Figure 4.14. In the case shown, s 
is increased uniformly from 0.2 to 1.2. Notice the substantial difference between the optimal performance 
and the optimal robust designs for large values of s. Other parametric curves can be easily obtained from 
Equations (4.17) and (4.18). 

4.1.6. Remarks. 

• The optimization of systems with uncertainty is a problem by itself. Developments and solutions 
based on the mere extension of existing deterministic tools tend to overlook the very primary nature 
of the problem. 

• In practice, each optimization problem has its own particular characteristics and objectives. A proper 
optimization problem must be posted by identifying the performance, robustness and reliability 
components of the desired outcome. Once the objective is set qualitatively, significant probabilistic 
indices must be designed. Only after these stages are consistently attained the optimization must 
be carried out. 

• The physical problem is deterministic in nature and the probabilistic approach is a mere consequence 
of the ignorance on the actual values of the uncertain parameters. 

• The Var[y] is proportional to the potential offset between the mean value of the output and the 
outcome of implementing the design. Then, variance reduction is a very important consideration in 
the study of engineering problems [9]. 

• Exact calculation of the moments of y (of any order) can be done if ft(x, a) is a polynomial function 
(of any order) of the random variables. If the dependence is not polynomial, a Taylor expansion or 
a interpolating function can be used to approximate n(x,a). Once such function is built, no further 
approximations are needed to derive explicit expressions for the moments of the output. Notice that 
nothing prevents us from using high order polynomials. In this framework, optimization problems 
can be easily solved. 

• If obtaining derivatives of Q(x, a) is difficult or computationally expensive, the sampling of few points 
on the interval of interest can be used to either (i) build an interpolating polynomial function or to 
(ii) approximate numerically the expressions for the moments. These practices were implemented 
and evaluated in Example (4.1.4). Problems with several inputs and several design variables can be 
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approached in the same manner [17], 

• The approximation of Q(x,a) generates an off-set between the true and the computed expressions 
for E[y ]. This error builds up as higher order moments of the output are calculated. 

• The reliability problem is the most intensive of all three due to the need for a detailed probabilistic 
description of the output. In the examples, the solution to this problem has been found by calculating 
fr(y) exactly and then evaluating the corresponding index. Obtaining exact and approximate 
expressions for the the density function of the output is in general very difficult. On the other 
hand, the use of a moments-based index seems to be insufficient to properly describe the output. 
Numerical methods such as FORM and SORM have been extensively used to study this problem 
[ 12 ]. 

• The reader must notice that the reliability problem has been intentionally precluded from Exam- 
ples (4.1.4) and (4.1.5). In general, information on the first two moments of the output is riot enough 
to properly build an accurate PDF. We could think of going further on our derivations and obtain 
the exact expressions for higher order moments and then use them to build a better PDF. While 
this practice leads to better approximations, an arbitrary distribution requires an infinite number of 
moments to be fully specified. This fact prevent us from studying the reliability problem confidently. 
Only if (i) the shape of the resultant PDF is known in advance or if (ii) the inclusion of additional 
moments in the construction of the PDF is inconsequent, truncation schemes such as the cumulant 
neglect closure method [6] should be used. 

• The approach used in Example (4.1.5) is extendable to cases with several independent random 
variables, not necessarily Gaussian. For well-known PDFs, the corresponding moment generating 
function and the property of independence can be used to systematically apply the expected value 
operator to polynomial expressions. In this fashion, explicit functions for the moments ( E[y ] and 
Var [y] among them) as a function of the design variable and the parameters of the PDFs of the 
inputs can be derived. Therefore, combinations of uniform, exponential, Weibull, Beta and gamma 
input random variables (among others) can be considered at no additional expense. 

• The minimum variance problem tends to eliminate the dependence of the output on the changes of 
the random variables. If the dependence on the design parameter does not define a tradeoff among 
two or more input random variables, one might be tempted to infer the solution. Notice however, 
that even in very simple problems (as in Example (4.1.2)) the robust solution is by no means trivial. 

4.2. Optimization of Control Systems. In this section we study the stability of some linear time- 
invariant control systems whose closed loop transfer function has one or two poles. We assume that the 
control structure is set but not the value of its gain, that for our purposes is the design variable. 

Assume that the closed loop response of a dynamic system is given by: 

n 

^nj(z,a)t(; ( ' l) ' = 0, (4.19) 

t=0 

where a is a random variable distributed according to /^(a), w is the system response, w (n >' is the n th 
derivative of w with respect to time t and x 6 * is the design variable (gain). The response, given by 
w = e xt , has the following characteristic equation: 

n 

Y yfli(x,a) = 0 (4.20) 

1=0 
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Denote the zeros of this Equation as Ai, A 2 . ■ ■ A„ and define the output as: 


y = maz{ Re{ Ai }, Re{A 2 }, . . . , Re{A n }}, (4-21) 

This quantity, which is a function of the design variable and the uncertain parameter, will be used to measure 
the quality of the controlled response in the presence of uncertainty. By moving the PDFs of y to the negative 
portion of the y- axis the expected value of the time needed to stabilize the system is reduced. On the other 
hand, the portion of fy{ y) on the positive part of the y-axis is a measure of the likelihood of instability in 
the response (probability of failure). In this framework, we are interested in finding which are the control 
designs, i.e. gains, that are optimal in terms of performance, robustness and reliability. 

This formulation allows us to perform the stability analysis of dynamic systems without dealing with 
the time dependence of the response nor with the probabilistic description of all eigenvalues. In this context, 
we will derive the exact expression for f Y (y) and use it to determine the optimal gains. 

Notice that the stability analysis of the first order system w - U(x,a)w = 0 leads to y = A = ft(x,a). 
This problem has the same form as the ones studied in Section 4.1. Numerical studies on the stability of 
control systems with uncertainty can be found in [15, 14]. 

4.2.1. Linear Dynamic System with Real Eigenvalues. Consider the dynamical system: 

w + (x 2 + a 2 - a - x)w + (a 2 x 2 + ax - x 3 - a 3 )w = 0, (4.22) 

where a is a random variable with PDF f A (a) and x is the gain to be determined. The corresponding 
eigenvalues are given by Ai = x - a 2 and A 2 = a- x 2 . Notice that the eigenvalues are real for all possible 
values of a and x. The Routh stability test of the corresponding deterministic system leads to the stable 
ranges {|x| > a} n {x < a 2 } for a > 0 and {x < a 2 } otherwise. Notice that, the system’s response is not only 
quite sensitive to changes in the uncertain parameter but also might run into unstable regimes. The PDF 
of the output is given by: 

f Y (y) = C[s(y) +t{y) + u(y)] (4.23) 


s{y) = 


a 3 1)v />'( a 3) if v € [ 2 / 3 , 2 / 2 ] 
0 otherwise 


(4.24) 


t(y) = 


/. 4 ( 2 / +x 2 ) if y € [2/1,00] U [-00, y 3 ] 
0 otherwise 


(4.25) 


u{y) = 


( 


f A (y + x 2 ) - a[ 1)v / aM + a ( 3 1)y fA(a 3 ) if 2 / € [ 2 / 2 , 20 ] and max{oi,a 2 } € [0, oc] 
0 otherwise 


(4.26) 


where C is a normalization constant, 2/1 = Ai(0), 2/2 = max{A 2 (a 2 ), A 2 (ai)}, 2/3 = rran{A 2 (a 2 ), A 2 (oi)}, 
aj = — 1 — x, a 2 = x, «3 = — \Jx - y, = —a 3 , <23 - 1 /( 204 ) and 04 = —03 

A particular case with a -4 N{ 1.1, 0.6) and x £ [-1.5, 1.5] is presented next. The Routh stabil- 
ity test of the corresponding deterministic system, i.e. a = m = 1.1, leads to the stable range x € 
[_1 .5,-1 .04] U [1.04, 1.21]. The spectrum of PDFs of y is shown in Figure 4.15. Equation (4.23) can also 
be applied to calculate numerically the three optimization indices. Figure 4.16 shows how they vary with x. 


15 
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Fig. 4.15. Spectrum of PDFs of the output for a Gaussian input random variable 



Fig. 4.16. First two order moments of the output and probability of failure as functions of the gain 


FYom this information optimal designs are then extracted. In the calculations, yu m was set as the threshold 
for stability, i.e. yi, rn = 0. 

Notice the difference in the admissible designs provided by the Routh test and the performance based 
criterion. If a deterministic solution within x € [1.04, 1.21] is implemented, the system has a high likelihood 
of being unstable. Notice also the sharp differences in the variance and probability of failure (instability) for 
the designs based merely on the expected value of the response. It is evident the need for better probabilistic 
descriptions of the response in order to select the best design among all the ones that are deterministically 
admissible. 

In this particular example the best gain is x = —1.5, for which both performance and reliability are 
optimal. The best robust design, in spite of being mathematically coherent has no practical value due to the 
behavior of E[y{x* ro )] and P{y{x' ro ) > y lim }. 

4.2.2. Linear Dynamic System with Complex Eigenvalues. Consider the dynamical system: 

yj + (x 2 - a)w + (a 2 + b — x)w = 0, (4-27) 
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Fig. 4.17. Root locus of the system for different values of the uncertain parameter 


where a is a random variable with PDF f A (a), b is a known parameter and x is the gain to be determined. 
The corresponding characteristic equation is: 

Ai ,2 = (a - x 2 ± V15)/2, 

where the discriminant D, is given by D = x 4 — 3a 2 — 2 ax 2 4 4x — 46. The eigenvalues of this characteristic 
equation are in general complex numbers. Root locus of the poles for multiple values of a are shown in 
Figure 4.17. Notice the excursion of the response to instability and the high sensitivity to changes in the 
uncertain parameter. 

According to the Routh stability test with a taking the value of its expected value, the system is stable 
on {|x| > v^} H {x < a 2 4 6} for a > 0, and on x < a 2 4 6 otherwise. Notice that in this non-probabilistic 
stability analysis there is no way to discriminate among all the globally stable control designs. The PDF of 
the output is given by: 

fy(y) = C[s(j/) + t{y) + u{y)\ (4-28) 

3{y) = I ° 3 )v fA(a 3 ) if V e [j/3,1/2] (4 . 29) 

I 0 otherwise 


t(y) = 


2f A {2y + x 2 ) if y € [ 2 / 1 , 00 ] U [- 00 , y 3 ] 
0 otherwise 


(4.30) 


u{y) = 


2f A (2y + x 2 ) - a[ 1)y f A {a. 1 ) + 4 1)v f A (a 3 ) if y € [2/2, yi] 
0 otherwise 


(4.31) 


where C is a normalization constant, y\ — Ai(a), x }2 = Ai(ei 2 )> 2/3 — Ai(ai), a — (—x 2 4 \A)A v — x A +Zx 

3 6 , ai = (—x 2 4 2y/u) /3, «2 = (— x 2 — 2y/v)j%, 03 = (2/ + \/«)/2, o>4 — (2/ — \/^)/2> K — 4(x ~ x2 V ~ — 3y , 

a< 1} » =(k 4 2x 2 4 Zy)/2K and a[ l)v =(k- 2x 2 - 3y)/2*. Notice that 6 < - </ 3/4 guarantees the existence 

of complex eigenvalues for all possible values of x. 

A numerical example with a -4 JV(1.55,0.5), b = -0.95 and x € [-1.5, 1.5] is presented next. The 
stability analysis of the corresponding deterministic system leads to the stable range x G {[-1.5, -1.24] U 
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Fig. 4.18. First two order moments of the output and probability of failure as functions of the gain 


[1.24, 1.45]}. Using Equation (4.28), the spectrum of PDFs as well as the solution to the three optimization 
problems is obtained. As before, yu m is taken as the stability threshold. 

The dependence of the design variable on the three indices is shown in Figure 4.18. Notice that some of 
the proper deterministic designs have high probability of failure as well as positive expected values. From this 
Figure is clear that x = — 1.5 leads to the best response in terms of performance and reliability. Nonetheless, 
the response is not exempt from running into unstability. While in both examples such designs coincide this 
is in general not the case. 

4.2.3. Remarks. 

• Stability control is an application in which considerations on both performance and reliability are 
crucial. Fast decaying transient responses are desirable in practice. 

• The analysis presented provides an exhaustive description of the implications of a particular design on 
the system response. This permits to attain optimal designs with the desired practical performances. 

• From the stability point of view, the reliability based design with yu m = 0 provides the gain with 
the least probability of making the system unstable. This class of problems, besides having a clear 
physical interpretation leads to results with significant practical value. 

• In the presence of uncertainty a designer is forced, in effect, to take a gamble. Under such cir- 
cumstances, rather than naively hoping for the best or over-conservatively focusing on the worst, 
the right decision consist in the best possible design whether favorable or unfavorable operating 
conditions occur. In this example, negative values of E[y\ do not necessarily stabilize the system. 
In practice, a and y are deterministic quantities that might be away from their corresponding mean 
values. It is possible to have the fortune of stabilizing a system with a design with E[y] > 0 as well 
as having the misfortune of destabilizing it using with a gain with E[y] < 0. Unless P{y > 0} = 0, 
we can not claim with certainty that the response is stable. 

• Results can be validated by comparing the behavior of E[y] for infinitely concentrated input variables 
with the results of applying the Routh stability test to the corresponding deterministic system. For 
a single input random variable this can be stated as follows. 

If 5(x) = E[y] when a — > iV(m,s) such that s -C 1 and x = {:r E x\e{x) < 0} then x must coincide 
with the results of the Routh stability test when f li(x,a) = D t (x,m) V?‘ in Equation (4.19). In both 
control examples, this practice led to consistent results. 
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5. Conclusions. This paper studies the problem of optimization of systems with uncertainty for per- 
formance, robustness and reliability based designs. By solving problems that admit a closed form solution 
for the PDF of the output, the effects of the design variable are fully evaluated in the probabilistic sense. 

Applications to static optimization show that these criteria not only require different formulations but 
also lead to different designs. The existence of bounds in the PDFs of the output, the relation between the 
probabilistic optimal solution and its deterministic counterpart, the effects of having several extreina. the 
extension of existing tools to problems that do not admit a complete probabilistic description of the output 
and strategies to manage multiple uncertain parameters are some of the aspects explored. 

Applications to stability control show the value of incorporating uncertainty in the early stages of the 
design. In the examples it is shown how T purely deterministic tools that disregard the uncertainty might 
lead to erroneous designs whose consequences can not only be far from optimal but also catastrophic. For 
example, we show that deterministic stability tests based on the expected value of the input parameter 
provide results with high probability of running into unstable regimes. 

It is very important to describe the objectives in terms of a consistent probabilistic formulation. Every 
problem has its own particular features and objectives. Failing in describing the desired response in terms 
of meaningful statistical indices will lead to faulty designs. A proper formulation as well as a complete 
probabilistic description of the output /response are crucial to attain designs with practical value. 
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